TAI profiles

From PES, one can obtain the TAI profiles through the mytTAI::TAI function. The myTAI package, to my knowledge, does not enable the flexibility to plot the separate sexes, for example in one plot. Therefore, I will not use the myTAI::PlotSignature function.

Here, I show the TAI profiles for the developmental stages. For tissues of the adult (matSP) Fucus species, one can find them in 2.1 TAI tissues in Fucus matSP.

library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.2     ✔ readr     2.1.4
## ✔ forcats   1.0.0     ✔ stringr   1.5.0
## ✔ ggplot2   3.4.3     ✔ tibble    3.2.1
## ✔ lubridate 1.9.2     ✔ tidyr     1.3.0
## ✔ purrr     1.0.2     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(myTAI)

Load PES

Non-transformed

Fd_PES <-
  readr::read_csv(file = "data/Fd_PES.csv")
## Rows: 7898 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (8): Phylostratum, gamete, E1, E2, E3, E4, E5, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fd_PES_F <-
  readr::read_csv(file = "data/Fd_PES_F.csv")
## Rows: 7898 Columns: 3
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (2): Phylostratum, gamete
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fd_PES_M <-
  readr::read_csv(file = "data/Fd_PES_M.csv")
## Rows: 7898 Columns: 3
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (2): Phylostratum, gamete
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES <-
  readr::read_csv(file = "data/Fs_PES.csv")
## Rows: 8291 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (8): Phylostratum, gamete, 24H, 48H, 1w, 3w, 4w, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES_F <-
  readr::read_csv(file = "data/Fs_PES_F.csv")
## Rows: 8291 Columns: 4
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (3): Phylostratum, gamete, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES_M <-
  readr::read_csv(file = "data/Fs_PES_M.csv")
## Rows: 8291 Columns: 4
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (3): Phylostratum, gamete, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Ec_PES_32m <-
  readr::read_csv(file = "data/Ec_PES_32m.csv")
## Rows: 11571 Columns: 11
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr  (1): GeneID
## dbl (10): Phylostratum, meiospore, immGA, matGA, oldGA, gamete, earlyPSP, im...
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Ec_PES_25f <-
  readr::read_csv(file = "data/Ec_PES_25f.csv")
## Rows: 11571 Columns: 11
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr  (1): GeneID
## dbl (10): Phylostratum, meiospore, immGA, matGA, oldGA, gamete, earlyPSP, im...
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

sqrt-tranformed

Fd_PES.sqrt <-
  readr::read_csv(file = "data/Fd_PES.sqrt.csv")
## Rows: 7898 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (8): Phylostratum, gamete, E1, E2, E3, E4, E5, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fd_PES_F.sqrt <-
  readr::read_csv(file = "data/Fd_PES_F.sqrt.csv")
## Rows: 7898 Columns: 3
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (2): Phylostratum, V1
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fd_PES_M.sqrt <-
  readr::read_csv(file = "data/Fd_PES_M.sqrt.csv")
## Rows: 7898 Columns: 3
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (2): Phylostratum, V1
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES.sqrt <-
  readr::read_csv(file = "data/Fs_PES.sqrt.csv")
## Rows: 8291 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (8): Phylostratum, gamete, 24H, 48H, 1w, 3w, 4w, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES_F.sqrt <-
  readr::read_csv(file = "data/Fs_PES_F.sqrt.csv")
## Rows: 8291 Columns: 4
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (3): Phylostratum, gamete, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES_M.sqrt <-
  readr::read_csv(file = "data/Fs_PES_M.sqrt.csv")
## Rows: 8291 Columns: 4
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (3): Phylostratum, gamete, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Ec_PES_32m.sqrt <-
  readr::read_csv(file = "data/Ec_PES_32m.sqrt.csv")
## Rows: 11571 Columns: 11
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr  (1): GeneID
## dbl (10): Phylostratum, meiospore, immGA, matGA, oldGA, gamete, earlyPSP, im...
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Ec_PES_25f.sqrt <-
  readr::read_csv(file = "data/Ec_PES_25f.sqrt.csv")
## Rows: 11571 Columns: 11
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr  (1): GeneID
## dbl (10): Phylostratum, meiospore, immGA, matGA, oldGA, gamete, earlyPSP, im...
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

log2-tranformed

Fd_PES.log2 <-
  readr::read_csv(file = "data/Fd_PES.log2.csv")
## Rows: 7898 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (8): Phylostratum, gamete, E1, E2, E3, E4, E5, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fd_PES_F.log2 <-
  readr::read_csv(file = "data/Fd_PES_F.log2.csv")
## Rows: 7898 Columns: 3
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (2): Phylostratum, V1
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fd_PES_M.log2 <-
  readr::read_csv(file = "data/Fd_PES_M.log2.csv")
## Rows: 7898 Columns: 3
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (2): Phylostratum, V1
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES.log2 <-
  readr::read_csv(file = "data/Fs_PES.log2.csv")
## Rows: 8291 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (8): Phylostratum, gamete, 24H, 48H, 1w, 3w, 4w, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES_F.log2 <-
  readr::read_csv(file = "data/Fs_PES_F.log2.csv")
## Rows: 8291 Columns: 4
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (3): Phylostratum, gamete, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES_M.log2 <-
  readr::read_csv(file = "data/Fs_PES_M.log2.csv")
## Rows: 8291 Columns: 4
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (3): Phylostratum, gamete, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Ec_PES_32m.log2 <-
  readr::read_csv(file = "data/Ec_PES_32m.log2.csv")
## Rows: 11571 Columns: 11
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr  (1): GeneID
## dbl (10): Phylostratum, meiospore, immGA, matGA, oldGA, gamete, earlyPSP, im...
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Ec_PES_25f.log2 <-
  readr::read_csv(file = "data/Ec_PES_25f.log2.csv")
## Rows: 11571 Columns: 11
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr  (1): GeneID
## dbl (10): Phylostratum, meiospore, immGA, matGA, oldGA, gamete, earlyPSP, im...
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

rank-tranformed

Fd_PES.rank <-
  readr::read_csv(file = "data/Fd_PES.rank.csv")
## Rows: 7898 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (8): Phylostratum, gamete, E1, E2, E3, E4, E5, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fd_PES_F.rank <-
  readr::read_csv(file = "data/Fd_PES_F.rank.csv")
## Rows: 7898 Columns: 3
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (2): Phylostratum, V1
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fd_PES_M.rank <-
  readr::read_csv(file = "data/Fd_PES_M.rank.csv")
## Rows: 7898 Columns: 3
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (2): Phylostratum, V1
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES.rank <-
  readr::read_csv(file = "data/Fs_PES.rank.csv")
## Rows: 8291 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (8): Phylostratum, gamete, 24H, 48H, 1w, 3w, 4w, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES_F.rank <-
  readr::read_csv(file = "data/Fs_PES_F.rank.csv")
## Rows: 8291 Columns: 4
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (3): Phylostratum, gamete, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES_M.rank <-
  readr::read_csv(file = "data/Fs_PES_M.rank.csv")
## Rows: 8291 Columns: 4
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (3): Phylostratum, gamete, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Ec_PES_32m.rank <-
  readr::read_csv(file = "data/Ec_PES_32m.rank.csv")
## Rows: 11571 Columns: 11
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr  (1): GeneID
## dbl (10): Phylostratum, meiospore, immGA, matGA, oldGA, gamete, earlyPSP, im...
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Ec_PES_25f.rank <-
  readr::read_csv(file = "data/Ec_PES_25f.rank.csv")
## Rows: 11571 Columns: 11
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr  (1): GeneID
## dbl (10): Phylostratum, meiospore, immGA, matGA, oldGA, gamete, earlyPSP, im...
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

rlog-tranformed

Fd_PES.rlog <-
  readr::read_csv(file = "data/Fd_PES.rlog.csv")
## Rows: 7898 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (8): Phylostratum, gamete, E1, E2, E3, E4, E5, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fd_PES_F.rlog <-
  readr::read_csv(file = "data/Fd_PES_F.rlog.csv")
## Rows: 7898 Columns: 3
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (2): Phylostratum, gamete
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fd_PES_M.rlog <-
  readr::read_csv(file = "data/Fd_PES_M.rlog.csv")
## Rows: 7898 Columns: 3
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (2): Phylostratum, gamete
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES.rlog <-
  readr::read_csv(file = "data/Fs_PES.rlog.csv")
## Rows: 8291 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (8): Phylostratum, gamete, 24H, 48H, 1w, 3w, 4w, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES_F.rlog <-
  readr::read_csv(file = "data/Fs_PES_F.rlog.csv")
## Rows: 8291 Columns: 4
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (3): Phylostratum, gamete, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES_M.rlog <-
  readr::read_csv(file = "data/Fs_PES_M.rlog.csv")
## Rows: 8291 Columns: 4
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (3): Phylostratum, gamete, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Ec_PES_32m.rlog <-
  readr::read_csv(file = "data/Ec_PES_32m.rlog.csv")
## Rows: 11571 Columns: 11
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr  (1): GeneID
## dbl (10): Phylostratum, meiospore, immGA, matGA, oldGA, gamete, earlyPSP, im...
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Ec_PES_25f.rlog <-
  readr::read_csv(file = "data/Ec_PES_25f.rlog.csv")
## Rows: 11571 Columns: 11
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr  (1): GeneID
## dbl (10): Phylostratum, meiospore, immGA, matGA, oldGA, gamete, earlyPSP, im...
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

myTAI issue

I will post an issue on GitHub but myTAI::tf() loses the column name hwen transforming only one column.

TAI calculation

I had previously used a function to combine males and females into one TAI plot.

#TI stands for transcriptome index
TI.preplot <- function(ExpressionSet, permutations = 50000){
  std <- 
    myTAI::bootMatrix(
      ExpressionSet = ExpressionSet, 
      permutations  = permutations) %>% 
    apply(2, stats::sd)
  TI.out <- 
    myTAI::TAI(ExpressionSet) %>%
    tibble::as_tibble(rownames = "Stage", colnames = c("TAI", "PS")) %>% 
    dplyr::bind_cols(as_tibble(std)) %>%
    dplyr::rename(TAI = 2, sd = 3)
  return(TI.out)
}

Fucus distichus

# function to process Fd. This will be changed to accommodate more seq data.
get_TAI_Fd <- function(PES_all, PES_M, PES_F, ordered_stages){
  TAI_b <- 
    TI.preplot(
      dplyr::select(PES_all, !c("gamete"))) %>%
    dplyr::mutate(Sex = "Mixed")
  stages_to_NA <- # this will differ for Fucus serratus
    ordered_stages[-c(1, 1+1)]
  TAI_M <- 
    TI.preplot(
      PES_M) %>%
    tibble::add_row(
      dplyr::filter(
        dplyr::select(TAI_b, !Sex), 
        Stage == ordered_stages[1+1])) %>%
    tibble::add_row(
      tibble::tibble(Stage = stages_to_NA, TAI = NA, sd = NA)
    ) %>%
    dplyr::mutate(Sex = "Male")
  TAI_F <- 
    TI.preplot(
      PES_F) %>%
    tibble::add_row(
      dplyr::filter(
        dplyr::select(TAI_b, !Sex), 
        Stage == ordered_stages[1+1])) %>%
    tibble::add_row(
      tibble::tibble(Stage = stages_to_NA, TAI = NA, sd = NA)
    ) %>%
    dplyr::mutate(Sex = "Female") # This is not true - this is needed for the plotting.
  TAI_out <- 
    dplyr::bind_rows(TAI_b, TAI_M, TAI_F)
  TAI_out$Stage <- base::factor(TAI_out$Stage, ordered_stages)
  return(TAI_out)
}
ordered_stages <- colnames(Fd_PES)[3:ncol(Fd_PES)]
Fd_TAI <- 
  get_TAI_Fd(
    PES_all = Fd_PES, 
    PES_M = Fd_PES_M, 
    PES_F = Fd_PES_F, 
    ordered_stages)
## New names:
## New names:
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`
Fd_TAI.log2 <- 
  get_TAI_Fd(
    PES_all = Fd_PES.log2, 
    PES_M = dplyr::rename(Fd_PES_M.log2, gamete = V1), 
    PES_F = dplyr::rename(Fd_PES_F.log2, gamete = V1), 
    ordered_stages)
## New names:
## New names:
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`
Fd_TAI.sqrt <- 
  get_TAI_Fd(
    PES_all = Fd_PES.sqrt, 
    PES_M = dplyr::rename(Fd_PES_M.sqrt, gamete = V1), 
    PES_F = dplyr::rename(Fd_PES_F.sqrt, gamete = V1), 
    ordered_stages)
## New names:
## New names:
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`
Fd_TAI.rank <- 
  get_TAI_Fd(
    PES_all = Fd_PES.rank, 
    PES_M = dplyr::rename(Fd_PES_M.rank, gamete = V1), 
    PES_F = dplyr::rename(Fd_PES_F.rank, gamete = V1), 
    ordered_stages)
## New names:
## New names:
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`
Fd_TAI.rlog <- 
  get_TAI_Fd(
    PES_all = Fd_PES.rlog, 
    PES_M = Fd_PES_M.rlog, 
    PES_F = Fd_PES_F.rlog, 
    ordered_stages)
## New names:
## New names:
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`

As I previously noted dplyr::rename() is used due to a bug I found in myTAI.

Fucus serratus

# function to process Fd. This will be changed to accommodate more seq data.
get_TAI_Fs <- function(PES_all, PES_M, PES_F, ordered_stages){
  TAI_b <- 
    TI.preplot(
      dplyr::select(PES_all, !c("gamete","matSP"))) %>%
    dplyr::mutate(Sex = "Mixed")
  stages_to_NA <- # this will differ for Fucus distichus
    ordered_stages[-c(1, 1+1, length(ordered_stages) - 1, length(ordered_stages))]
  TAI_M <- 
    TI.preplot(
      PES_M) %>%
    tibble::add_row(
      dplyr::filter(
        dplyr::select(TAI_b, !Sex), 
        Stage %in% c(
          ordered_stages[1+1], 
          ordered_stages[length(ordered_stages) - 1]))) %>%
    tibble::add_row(
      tibble::tibble(Stage = stages_to_NA, TAI = NA, sd = NA)
    ) %>%
    dplyr::mutate(Sex = "Male")
  TAI_F <- 
    TI.preplot(
      PES_F) %>%
    tibble::add_row(
      dplyr::filter(
        dplyr::select(TAI_b, !Sex), 
        Stage %in% c(
          ordered_stages[1+1], 
          ordered_stages[length(ordered_stages) - 1]))) %>%
    tibble::add_row(
      tibble::tibble(Stage = stages_to_NA, TAI = NA, sd = NA)
    ) %>%
    dplyr::mutate(Sex = "Female") # This is not true. this is needed for the plotting.
  TAI_out <- 
    dplyr::bind_rows(TAI_b, TAI_M, TAI_F)
  TAI_out$Stage <- base::factor(TAI_out$Stage, ordered_stages)
  return(TAI_out)
}
ordered_stages <- colnames(Fs_PES)[3:ncol(Fs_PES)]
Fs_TAI <- 
  get_TAI_Fs(
    PES_all = Fs_PES, 
    PES_M = Fs_PES_M, 
    PES_F = Fs_PES_F, 
    ordered_stages)
## New names:
## New names:
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`
Fs_TAI.log2 <- 
  get_TAI_Fs(
    PES_all = Fs_PES.log2, 
    PES_M = Fs_PES_M.log2, 
    PES_F = Fs_PES_F.log2,
    ordered_stages)
## New names:
## New names:
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`
Fs_TAI.sqrt <- 
  get_TAI_Fs(
    PES_all = Fs_PES.sqrt, 
    PES_M = Fs_PES_M.sqrt, 
    PES_F = Fs_PES_F.sqrt, 
    ordered_stages)
## New names:
## New names:
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`
Fs_TAI.rank <- 
  get_TAI_Fs(
    PES_all = Fs_PES.rank, 
    PES_M = Fs_PES_M.rank, 
    PES_F = Fs_PES_F.rank, 
    ordered_stages)
## New names:
## New names:
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`
Fs_TAI.rlog <- 
  get_TAI_Fs(
    PES_all = Fs_PES.rlog, 
    PES_M = Fs_PES_M.rlog, 
    PES_F = Fs_PES_F.rlog, 
    ordered_stages)
## New names:
## New names:
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`

As I previously noted dplyr::rename() is used due to a bug I found in myTAI.

Ectocarpus

# function to process Fd. This will be changed to accommodate more seq data.
get_TAI_Ec <- function(PES_M, PES_F, ordered_stages){
  TAI_M <- 
    TI.preplot(PES_M) %>%
    dplyr::mutate(Sex = "Male")
  TAI_F <- 
    TI.preplot(PES_F) %>%
    dplyr::mutate(Sex = "Female")
  stages_to_NA <- # this will differ for Fucus distichus
    ordered_stages[-c(1, 1+1, length(ordered_stages) - 1, length(ordered_stages))]
  TAI_out <- 
    dplyr::bind_rows(TAI_M, TAI_F)
  TAI_out$Stage <- base::factor(TAI_out$Stage, ordered_stages)
  return(TAI_out)
}
ordered_stages <- colnames(Ec_PES_25f)[3:ncol(Ec_PES_25f)]

Ec_TAI <- 
  get_TAI_Ec(
    PES_M = Ec_PES_32m, 
    PES_F = Ec_PES_25f, 
    ordered_stages)
## New names:
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`
Ec_TAI.sqrt <- 
  get_TAI_Ec(
    PES_M = Ec_PES_32m.sqrt, 
    PES_F = Ec_PES_25f.sqrt, 
    ordered_stages)
## New names:
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`
Ec_TAI.log2 <- 
  get_TAI_Ec(
    PES_M = Ec_PES_32m.log2, 
    PES_F = Ec_PES_25f.log2, 
    ordered_stages)
## New names:
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`
Ec_TAI.rank <- 
  get_TAI_Ec(
    PES_M = Ec_PES_32m.rank, 
    PES_F = Ec_PES_25f.rank, 
    ordered_stages)
## New names:
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`
Ec_TAI.rlog <- 
  get_TAI_Ec(
    PES_M = Ec_PES_32m.rlog, 
    PES_F = Ec_PES_25f.rlog, 
    ordered_stages)
## New names:
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`

TAI profiles

Fucus serratus

Fs_TAI %>% 
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = Sex,
      colour = Sex,
      fill = Sex)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), alpha = 0.1) +
  ggplot2::geom_line(size = 2, lineend = "round") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus serratus",
       subtitle = "TPM")
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

Fs_TAI.sqrt %>% 
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = Sex,
      colour = Sex,
      fill = Sex)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), alpha = 0.1) +
  ggplot2::geom_line(size = 2, lineend = "round") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus serratus",
       subtitle = "sqrt(TPM)")

Fs_TAI.log2 %>% 
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = Sex,
      colour = Sex,
      fill = Sex)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), alpha = 0.1) +
  ggplot2::geom_line(size = 2, lineend = "round") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus serratus",
       subtitle = "log2(TPM+\u03b1)")

Fs_TAI.rank %>% 
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = Sex,
      colour = Sex,
      fill = Sex)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), alpha = 0.1) +
  ggplot2::geom_line(size = 2, lineend = "round") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus serratus",
       subtitle = "rank(TPM)")

Fs_TAI.rlog %>% 
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = Sex,
      colour = Sex,
      fill = Sex)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), alpha = 0.1) +
  ggplot2::geom_line(size = 2, lineend = "round") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus serratus",
       subtitle = "rlog(TPM)")

The dynamic range of the rlog transform is really small. Does this mean that the signal isn’t biologically strong?

Fucus distichus

Fd_TAI %>% 
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = Sex,
      colour = Sex,
      fill = Sex)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), alpha = 0.1) +
  ggplot2::geom_line(size = 2, lineend = "round") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus distichus",
       subtitle = "TPM")
## Warning: Removed 10 rows containing missing values (`geom_line()`).

Fd_TAI.sqrt %>% 
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = Sex,
      colour = Sex,
      fill = Sex)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), alpha = 0.1) +
  ggplot2::geom_line(size = 2, lineend = "round") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus distichus",
       subtitle = "sqrt(TPM)")
## Warning: Removed 10 rows containing missing values (`geom_line()`).

Fd_TAI.log2 %>% 
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = Sex,
      colour = Sex,
      fill = Sex)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), alpha = 0.1) +
  ggplot2::geom_line(size = 2, lineend = "round") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus distichus",
       subtitle = "log2(TPM+\u03b1)")
## Warning: Removed 10 rows containing missing values (`geom_line()`).

Fd_TAI.rank %>% 
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = Sex,
      colour = Sex,
      fill = Sex)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), alpha = 0.1) +
  ggplot2::geom_line(size = 2, lineend = "round") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus distichus",
       subtitle = "rank(TPM)")
## Warning: Removed 10 rows containing missing values (`geom_line()`).

Fd_TAI.rlog %>% 
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = Sex,
      colour = Sex,
      fill = Sex)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), alpha = 0.1) +
  ggplot2::geom_line(size = 2, lineend = "round") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus distichus",
       subtitle = "rlog(TPM)")
## Warning: Removed 10 rows containing missing values (`geom_line()`).

TAI Fucus embryonic stages only

#Fs
Fs_TAI %>% 
  dplyr::filter(Stage %in% c("24H", "48H", "1w", "3w", "4w") & Sex == "Mixed") %>%
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = 1)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), 
    colour = "#00A087FF", 
    fill = "#00A087FF", 
    alpha = 0.1) +
  ggplot2::geom_line(
    size = 2,
    lineend = "round",
    colour = "#00A087FF") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus serratus",
       subtitle = "TPM")

Fs_TAI.sqrt %>% 
  dplyr::filter(Stage %in% c("24H", "48H", "1w", "3w", "4w") & Sex == "Mixed") %>%
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = 1)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), 
    colour = "#00A087FF", 
    fill = "#00A087FF", 
    alpha = 0.1) +
  ggplot2::geom_line(
    size = 2,
    lineend = "round",
    colour = "#00A087FF") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus serratus",
       subtitle = "sqrt(TPM)")

Fs_TAI.log2 %>% 
  dplyr::filter(Stage %in% c("24H", "48H", "1w", "3w", "4w") & Sex == "Mixed") %>%
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = 1)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), 
    colour = "#00A087FF", 
    fill = "#00A087FF", 
    alpha = 0.1) +
  ggplot2::geom_line(
    size = 2,
    lineend = "round",
    colour = "#00A087FF") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus serratus",
       subtitle = "log2(TPM+\u03b1)")

Fs_TAI.rank %>% 
  dplyr::filter(Stage %in% c("24H", "48H", "1w", "3w", "4w") & Sex == "Mixed") %>%
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = 1)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), 
    colour = "#00A087FF", 
    fill = "#00A087FF", 
    alpha = 0.1) +
  ggplot2::geom_line(
    size = 2,
    lineend = "round",
    colour = "#00A087FF") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus serratus",
       subtitle = "rank(TPM)")

Fs_TAI.rlog %>% 
  dplyr::filter(Stage %in% c("24H", "48H", "1w", "3w", "4w") & Sex == "Mixed") %>%
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = 1)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), 
    colour = "#00A087FF", 
    fill = "#00A087FF", 
    alpha = 0.1) +
  ggplot2::geom_line(
    size = 2,
    lineend = "round",
    colour = "#00A087FF") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus serratus",
       subtitle = "rlog(TPM)")

# Fd
Fd_TAI %>% 
  dplyr::filter(Stage %in% c("E1", "E2", "E3", "E4", "E5") & Sex == "Mixed") %>%
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = 1)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), 
    colour = "#00A087FF", 
    fill = "#00A087FF", 
    alpha = 0.1) +
  ggplot2::geom_line(
    size = 2,
    lineend = "round",
    colour = "#00A087FF") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus distichus",
       subtitle = "TPM")

Fd_TAI.sqrt %>% 
  dplyr::filter(Stage %in% c("E1", "E2", "E3", "E4", "E5") & Sex == "Mixed") %>%
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = 1)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), 
    colour = "#00A087FF", 
    fill = "#00A087FF", 
    alpha = 0.1) +
  ggplot2::geom_line(
    size = 2,
    lineend = "round",
    colour = "#00A087FF") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus distichus",
       subtitle = "sqrt(TPM)")

Fd_TAI.log2 %>% 
  dplyr::filter(Stage %in% c("E1", "E2", "E3", "E4", "E5") & Sex == "Mixed") %>%
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = 1)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), 
    colour = "#00A087FF", 
    fill = "#00A087FF", 
    alpha = 0.1) +
  ggplot2::geom_line(
    size = 2,
    lineend = "round",
    colour = "#00A087FF") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus distichus",
       subtitle = "log2(TPM+\u03b1)")

Fd_TAI.rank %>% 
  dplyr::filter(Stage %in% c("E1", "E2", "E3", "E4", "E5") & Sex == "Mixed") %>%
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = 1)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), 
    colour = "#00A087FF", 
    fill = "#00A087FF", 
    alpha = 0.1) +
  ggplot2::geom_line(
    size = 2,
    lineend = "round",
    colour = "#00A087FF") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus distichus",
       subtitle = "rank(TPM)")

Fd_TAI.rlog %>% 
  dplyr::filter(Stage %in% c("E1", "E2", "E3", "E4", "E5") & Sex == "Mixed") %>%
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = 1)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), 
    colour = "#00A087FF", 
    fill = "#00A087FF", 
    alpha = 0.1) +
  ggplot2::geom_line(
    size = 2,
    lineend = "round",
    colour = "#00A087FF") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus distichus",
       subtitle = "rlog(TPM)")

Ectocarpus

Ec_TAI %>% 
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = Sex,
      colour = Sex,
      fill = Sex)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), alpha = 0.1) +
  ggplot2::geom_line(size = 2, lineend = "round") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  ggplot2::scale_x_discrete(guide = ggplot2::guide_axis(n.dodge=2)) +
  labs(title = "Ectocarpus",
       subtitle = "TPM")

Ec_TAI.sqrt %>% 
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = Sex,
      colour = Sex,
      fill = Sex)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), alpha = 0.1) +
  ggplot2::geom_line(size = 2, lineend = "round") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  ggplot2::scale_x_discrete(guide = ggplot2::guide_axis(n.dodge=2)) +
  labs(title = "Ectocarpus",
       subtitle = "sqrt(TPM)")

Ec_TAI.log2 %>% 
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = Sex,
      colour = Sex,
      fill = Sex)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), alpha = 0.1) +
  ggplot2::geom_line(size = 2, lineend = "round") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  ggplot2::scale_x_discrete(guide = ggplot2::guide_axis(n.dodge=2)) +
  labs(title = "Ectocarpus",
       subtitle = "log2(TPM+\u03b1)")

Ec_TAI.rank %>% 
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = Sex,
      colour = Sex,
      fill = Sex)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), alpha = 0.1) +
  ggplot2::geom_line(size = 2, lineend = "round") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  ggplot2::scale_x_discrete(guide = ggplot2::guide_axis(n.dodge=2)) +
  labs(title = "Ectocarpus",
       subtitle = "rank(TPM)")

Ec_TAI.rlog %>% 
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = Sex,
      colour = Sex,
      fill = Sex)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), alpha = 0.1) +
  ggplot2::geom_line(size = 2, lineend = "round") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  ggplot2::scale_x_discrete(guide = ggplot2::guide_axis(n.dodge=2)) +
  labs(title = "Ectocarpus",
       subtitle = "rlog(TPM)")

The r-log is behaving weird because the transformation is non-linear but also does not preserve the monotonic relationship.

tfStability

Flat line test

I will later redo with 50 000 permutations but I will use 20 000 for now since I get an error Error in gamma_MME[[3]] : subscript out of bounds.

I have redone this and it works without error.

F. serratus

Fs_PES_tS <- 
  tfStability(
    dplyr::select(Fs_PES, !c(gamete, matSP)),
    transforms = c("none", "sqrt", "log2", "rank"),
    permutations = 50000
    )
Fs_PES_tS_processed <-
  tibble::as_tibble(Fs_PES_tS, rownames = "transformation")

Fs_PES_tS.rlog <- 
  tfStability(
    dplyr::select(Fs_PES.rlog, !c(gamete, matSP)),
    transforms = c("none"), # it is already transformed.
    permutations = 50000
    )
Fs_PES_tS.rlog_processed <-
  tibble::as_tibble(Fs_PES_tS.rlog, rownames = "transformation") %>%
  dplyr::mutate(transformation = "rlog")

Fs_PES_tS_res <- 
  dplyr::bind_rows(Fs_PES_tS_processed, Fs_PES_tS.rlog_processed) %>%
  dplyr::rename("pval" = value) %>%
  dplyr::mutate(test = "FlatLineTest")

F. distichus

Fd_PES_tS <- 
  tfStability(
    dplyr::select(Fd_PES, !c(gamete, matSP)),
    transforms = c("none", "sqrt", "log2", "rank"),
    permutations = 50000
    )
Fd_PES_tS_processed <-
  tibble::as_tibble(Fd_PES_tS, rownames = "transformation")

Fd_PES_tS.rlog <- 
  tfStability(
    dplyr::select(Fd_PES.rlog, !c(gamete, matSP)),
    transforms = c("none"), # it is already transformed.
    permutations = 50000
    )
Fd_PES_tS.rlog_processed <-
  tibble::as_tibble(Fd_PES_tS.rlog, rownames = "transformation") %>%
  dplyr::mutate(transformation = "rlog")

Fd_PES_tS_res <- 
  dplyr::bind_rows(Fd_PES_tS_processed, Fd_PES_tS.rlog_processed) %>%
  dplyr::rename("pval" = value) %>%
  dplyr::mutate(test = "FlatLineTest")

Ectocarpus

# female
Ec_PES_25f_tS <- 
  tfStability(
    Ec_PES_25f,
    transforms = c("none", "sqrt", "log2", "rank"),
    permutations = 50000
    )
Ec_PES_25f_tS_processed <-
  tibble::as_tibble(Ec_PES_25f_tS, rownames = "transformation")

Ec_PES_25f_tS.rlog <- 
  tfStability(
    Ec_PES_25f.rlog,
    transforms = c("none"), # it is already transformed.
    permutations = 50000
    )
Ec_PES_25f_tS.rlog_processed <-
  tibble::as_tibble(Ec_PES_25f_tS.rlog, rownames = "transformation") %>%
  dplyr::mutate(transformation = "rlog")

Ec_PES_25f_tS_res <- 
  dplyr::bind_rows(Ec_PES_25f_tS_processed, Ec_PES_25f_tS.rlog_processed) %>%
  dplyr::rename("pval" = value) %>%
  dplyr::mutate(test = "FlatLineTest")

# male
Ec_PES_32m_tS <- 
  tfStability(
    Ec_PES_32m,
    transforms = c("none", "sqrt", "log2", "rank"),
    permutations = 50000
    )
Ec_PES_32m_tS_processed <-
  tibble::as_tibble(Ec_PES_32m_tS, rownames = "transformation")

Ec_PES_32m_tS.rlog <- 
  tfStability(
    Ec_PES_32m.rlog,
    transforms = c("none"), # it is already transformed.
    permutations = 50000
    )
Ec_PES_32m_tS.rlog_processed <-
  tibble::as_tibble(Ec_PES_32m_tS.rlog, rownames = "transformation") %>%
  dplyr::mutate(transformation = "rlog")

Ec_PES_32m_tS_res <- 
  dplyr::bind_rows(Ec_PES_32m_tS_processed, Ec_PES_32m_tS.rlog_processed) %>%
  dplyr::rename("pval" = value) %>%
  dplyr::mutate(test = "FlatLineTest")

Plot p-values

Fucus serratus development

Fs_PES_tS_res %>%
  ggplot2::ggplot(
    aes(
      y = -log10(pval), 
      x = transformation)) +
  ggplot2::geom_col(width = 0.05) +
  ggplot2::geom_point(size = 5) +
  ggplot2::geom_hline(
    yintercept = -log10(0.05), 
    colour = "#E64B35FF",
    linetype = 'dashed',
    size = 2
    ) +
  ggplot2::annotate(
    "label",
    x = 0.7, 
    y = max(-log10(Fs_PES_tS_res$pval))/6,
    label = "pval < 0.05",
    # angle = 270,
    colour = "#E64B35FF") +
  ggplot2::labs(
    title = "Fucus serratus",
    subtitle = "FlatLineTest",
    x = "Transformation"
  ) +
  ggplot2::coord_flip() +
  ggplot2::theme_classic()

# or perhaps we can use raw Pvalues
Fs_PES_tS_res %>%
  ggplot2::ggplot(
    aes(
      y = pval, 
      x = transformation)) +
  ggplot2::geom_col(width = 0.05) +
  ggplot2::geom_point(size = 5) +
  ggplot2::geom_hline(
    yintercept = 0.05, 
    colour = "#E64B35FF",
    linetype = 'dashed',
    size = 2
    ) +
  ggplot2::annotate(
    "label",
    x = 1.5, 
    y = 0.05 - 0.005,
    label = "pval < 0.05",
    # angle = 270,
    colour = "#E64B35FF") +
  ggplot2::labs(
    title = "Fucus serratus",
    subtitle = "FlatLineTest",
    x = "Transformation"
  ) +
  ggplot2::coord_flip() +
  ggplot2::theme_classic()

Fucus distichus development

Fd_PES_tS_res %>%
  ggplot2::ggplot(
    aes(
      y = -log10(pval), 
      x = transformation)) +
  ggplot2::geom_col(width = 0.05) +
  ggplot2::geom_point(size = 5) +
  ggplot2::geom_hline(
    yintercept = -log10(0.05), 
    colour = "#E64B35FF",
    linetype = 'dashed',
    size = 2
    ) +
  ggplot2::annotate(
    "label",
    x = 0.7, 
    y = max(-log10(Fd_PES_tS_res$pval))/6,
    label = "pval < 0.05",
    # angle = 270,
    colour = "#E64B35FF") +
  ggplot2::labs(
    title = "Fucus distichus",
    subtitle = "FlatLineTest",
    x = "Transformation"
  ) +
  ggplot2::coord_flip() +
  ggplot2::theme_classic()

# or perhaps we can use raw Pvalues
Fd_PES_tS_res %>%
  ggplot2::ggplot(
    aes(
      y = pval, 
      x = transformation)) +
  ggplot2::geom_col(width = 0.05) +
  ggplot2::geom_point(size = 5) +
  ggplot2::geom_hline(
    yintercept = 0.05, 
    colour = "#E64B35FF",
    linetype = 'dashed',
    size = 2
    ) +
  ggplot2::annotate(
    "label",
    x = 1.5, 
    y = 0.05 - 0.005,
    label = "pval < 0.05",
    # angle = 270,
    colour = "#E64B35FF") +
  ggplot2::labs(
    title = "Fucus distichus",
    subtitle = "FlatLineTest",
    x = "Transformation"
  ) +
  ggplot2::coord_flip() +
  ggplot2::theme_classic()

Ec_PES_32m_tS_res %>%
  ggplot2::ggplot(
    aes(
      y = -log10(pval), 
      x = transformation)) +
  ggplot2::geom_col(width = 0.05) +
  ggplot2::geom_point(size = 5) +
  ggplot2::geom_hline(
    yintercept = -log10(0.05), 
    colour = "#E64B35FF",
    linetype = 'dashed',
    size = 2
    ) +
  ggplot2::annotate(
    "label",
    x = 0.7, 
    y = max(-log10(Ec_PES_32m_tS_res$pval))/6,
    label = "pval < 0.05",
    # angle = 270,
    colour = "#E64B35FF") +
  ggplot2::labs(
    title = "Ectocarpus (male)",
    subtitle = "FlatLineTest",
    x = "Transformation"
  ) +
  ggplot2::coord_flip() +
  ggplot2::theme_classic()

Ec_PES_32m_tS_res %>%
  ggplot2::ggplot(
    aes(
      y = pval, 
      x = transformation)) +
  ggplot2::geom_col(width = 0.05) +
  ggplot2::geom_point(size = 5) +
  ggplot2::geom_hline(
    yintercept = 0.05, 
    colour = "#E64B35FF",
    linetype = 'dashed',
    size = 2
    ) +
  ggplot2::annotate(
    "label",
    x = 1.5, 
    y = 0.05 - 0.005,
    label = "pval < 0.05",
    # angle = 270,
    colour = "#E64B35FF") +
  ggplot2::labs(
    title = "Ectocarpus (male)",
    subtitle = "FlatLineTest",
    x = "Transformation"
  ) +
  ggplot2::coord_flip() +
  ggplot2::theme_classic()

Ec_PES_25f_tS_res %>%
  ggplot2::ggplot(
    aes(
      y = -log10(pval), 
      x = transformation)) +
  ggplot2::geom_col(width = 0.05) +
  ggplot2::geom_point(size = 5) +
  ggplot2::geom_hline(
    yintercept = -log10(0.05), 
    colour = "#E64B35FF",
    linetype = 'dashed',
    size = 2
    ) +
  ggplot2::annotate(
    "label",
    x = 0.7, 
    y = max(-log10(Ec_PES_25f_tS_res$pval))/6,
    label = "pval < 0.05",
    # angle = 270,
    colour = "#E64B35FF") +
  ggplot2::labs(
    title = "Ectocarpus (female)",
    subtitle = "FlatLineTest",
    x = "Transformation"
  ) +
  ggplot2::coord_flip() +
  ggplot2::theme_classic()

Ec_PES_25f_tS_res %>%
  ggplot2::ggplot(
    aes(
      y = pval, 
      x = transformation)) +
  ggplot2::geom_col(width = 0.05) +
  ggplot2::geom_point(size = 5) +
  ggplot2::geom_hline(
    yintercept = 0.05, 
    colour = "#E64B35FF",
    linetype = 'dashed',
    size = 2
    ) +
  ggplot2::annotate(
    "label",
    x = 1.5, 
    y = 0.05 - 0.005,
    label = "pval < 0.05",
    # angle = 270,
    colour = "#E64B35FF") +
  ggplot2::labs(
    title = "Ectocarpus (female)",
    subtitle = "FlatLineTest",
    x = "Transformation"
  ) +
  ggplot2::coord_flip() +
  ggplot2::theme_classic()

Non-flat line profile

Fs_PES_tS <- 
  tfStability(
    dplyr::select(Fs_PES, !c(gamete, matSP)),
    TestStatistic = "ReductiveHourglassTest",
    transforms = c("none", "sqrt", "log2", "rank"),
    modules = list(early = 1:2, mid = 3:4, late = 5),
    permutations = 50000
    )
## Proceeding with the ReductiveHourglassTest
Fs_PES_tS_processed <-
  tibble::as_tibble(Fs_PES_tS, rownames = "transformation")

Fs_PES_tS.rlog <- 
  tfStability(
    dplyr::select(Fs_PES.rlog, !c(gamete, matSP)),
    TestStatistic = "ReductiveHourglassTest",
    transforms = c("none"), # it is already transformed.
    modules = list(early = 1:2, mid = 3:4, late = 5),
    permutations = 50000
    )
## Proceeding with the ReductiveHourglassTest
Fs_PES_tS.rlog_processed <-
  tibble::as_tibble(Fs_PES_tS.rlog, rownames = "transformation") %>%
  dplyr::mutate(transformation = "rlog")

Fs_PES_tS_res <- 
  dplyr::bind_rows(Fs_PES_tS_processed, Fs_PES_tS.rlog_processed) %>%
  dplyr::rename("pval" = value) %>%
  dplyr::mutate(test = "ReductiveHourglassTest")

Fd_PES_tS <- 
  tfStability(
    dplyr::select(Fd_PES, !c(gamete, matSP)),
    TestStatistic = "LateConservationTest",
    transforms = c("none", "sqrt", "log2", "rank"),
    modules = list(early = 1:2, mid = 3:4, late = 5),
    permutations = 50000
    )
## Proceeding with the LateConservationTest
Fd_PES_tS_processed <-
  tibble::as_tibble(Fd_PES_tS, rownames = "transformation")

Fd_PES_tS.rlog <- 
  tfStability(
    dplyr::select(Fd_PES.rlog, !c(gamete, matSP)),
    TestStatistic = "LateConservationTest",
    transforms = c("none"), # it is already transformed.
    modules = list(early = 1:2, mid = 3:4, late = 5),
    permutations = 50000
    )
## Proceeding with the LateConservationTest
Fd_PES_tS.rlog_processed <-
  tibble::as_tibble(Fd_PES_tS.rlog, rownames = "transformation") %>%
  dplyr::mutate(transformation = "rlog")

Fd_PES_tS_res <- 
  dplyr::bind_rows(Fd_PES_tS_processed, Fd_PES_tS.rlog_processed) %>%
  dplyr::rename("pval" = value) %>%
  dplyr::mutate(test = "LateConservationTest")

Plot p-values

Fucus serratus development

Fs_PES_tS_res %>%
  ggplot2::ggplot(
    aes(
      y = -log10(pval), 
      x = transformation)) +
  ggplot2::geom_col(width = 0.05) +
  ggplot2::geom_point(size = 5) +
  ggplot2::geom_hline(
    yintercept = -log10(0.05), 
    colour = "#E64B35FF",
    linetype = 'dashed',
    size = 2
    ) +
  ggplot2::annotate(
    "label",
    x = 0.7, 
    y = max(-log10(Fs_PES_tS_res$pval))/6,
    label = "pval < 0.05",
    # angle = 270,
    colour = "#E64B35FF") +
  ggplot2::labs(
    title = "Fucus serratus",
    subtitle = "ReductiveHourglassTest",
    x = "Transformation"
  ) +
  ggplot2::coord_flip() +
  ggplot2::theme_classic()

# or perhaps we can use raw Pvalues
Fs_PES_tS_res %>%
  ggplot2::ggplot(
    aes(
      y = pval, 
      x = transformation)) +
  ggplot2::geom_col(width = 0.05) +
  ggplot2::geom_point(size = 5) +
  ggplot2::geom_hline(
    yintercept = 0.05, 
    colour = "#E64B35FF",
    linetype = 'dashed',
    size = 2
    ) +
  ggplot2::annotate(
    "label",
    x = 1.5, 
    y = 0.05 - 0.005,
    label = "pval < 0.05",
    # angle = 270,
    colour = "#E64B35FF") +
  ggplot2::labs(
    title = "Fucus serratus",
    subtitle = "ReductiveHourglassTest",
    x = "Transformation"
  ) +
  ggplot2::coord_flip() +
  ggplot2::theme_classic()

Fucus distichus development

Fd_PES_tS_res %>%
  ggplot2::ggplot(
    aes(
      y = -log10(pval), 
      x = transformation)) +
  ggplot2::geom_col(width = 0.05) +
  ggplot2::geom_point(size = 5) +
  ggplot2::geom_hline(
    yintercept = -log10(0.05), 
    colour = "#E64B35FF",
    linetype = 'dashed',
    size = 2
    ) +
  ggplot2::annotate(
    "label",
    x = 0.7, 
    y = max(-log10(Fd_PES_tS_res$pval))/6,
    label = "pval < 0.05",
    # angle = 270,
    colour = "#E64B35FF") +
  ggplot2::labs(
    title = "Fucus distichus",
    subtitle = "LateConservationTest",
    x = "Transformation"
  ) +
  ggplot2::coord_flip() +
  ggplot2::theme_classic()

# y = -log10(0.05) + max(Fd_PES_tS_res)/9,
# or perhaps we can use raw Pvalues
Fd_PES_tS_res %>%
  ggplot2::ggplot(
    aes(
      y = pval, 
      x = transformation)) +
  ggplot2::geom_col(width = 0.05) +
  ggplot2::geom_point(size = 5) +
  ggplot2::geom_hline(
    yintercept = 0.05, 
    colour = "#E64B35FF",
    linetype = 'dashed',
    size = 2
    ) +
  ggplot2::annotate(
    "label",
    x = 1.5, 
    y = 0.05 - 0.005,
    label = "pval < 0.05",
    # angle = 270,
    colour = "#E64B35FF") +
  ggplot2::labs(
    title = "Fucus distichus",
    subtitle = "LateConservationTest",
    x = "Transformation"
  ) +
  ggplot2::coord_flip() +
  ggplot2::theme_classic()

Get session info.

devtools::session_info()
## ─ Session info ───────────────────────────────────────────────────────────────
##  setting  value
##  version  R version 4.2.2 (2022-10-31)
##  os       macOS Big Sur ... 10.16
##  system   x86_64, darwin17.0
##  ui       X11
##  language (EN)
##  collate  en_GB.UTF-8
##  ctype    en_GB.UTF-8
##  tz       Europe/Berlin
##  date     2023-08-29
##  pandoc   2.19.2 @ /Applications/RStudio.app/Contents/MacOS/quarto/bin/tools/ (via rmarkdown)
## 
## ─ Packages ───────────────────────────────────────────────────────────────────
##  package      * version    date (UTC) lib source
##  bit            4.0.5      2022-11-15 [1] CRAN (R 4.2.0)
##  bit64          4.0.5      2020-08-30 [1] CRAN (R 4.2.0)
##  bslib          0.5.1      2023-08-11 [1] CRAN (R 4.2.0)
##  cachem         1.0.8      2023-05-01 [1] CRAN (R 4.2.0)
##  callr          3.7.3      2022-11-02 [1] CRAN (R 4.2.0)
##  cli            3.6.1      2023-03-23 [1] CRAN (R 4.2.0)
##  codetools      0.2-19     2023-02-01 [1] CRAN (R 4.2.0)
##  colorspace     2.1-0      2023-01-23 [1] CRAN (R 4.2.0)
##  crayon         1.5.2      2022-09-29 [1] CRAN (R 4.2.0)
##  devtools       2.4.5      2022-10-11 [1] CRAN (R 4.2.0)
##  digest         0.6.33     2023-07-07 [1] CRAN (R 4.2.0)
##  dplyr        * 1.1.2      2023-04-20 [1] CRAN (R 4.2.0)
##  ellipsis       0.3.2      2021-04-29 [1] CRAN (R 4.2.0)
##  evaluate       0.21       2023-05-05 [1] CRAN (R 4.2.0)
##  fansi          1.0.4      2023-01-22 [1] CRAN (R 4.2.0)
##  farver         2.1.1      2022-07-06 [1] CRAN (R 4.2.0)
##  fastmap        1.1.1      2023-02-24 [1] CRAN (R 4.2.0)
##  fitdistrplus   1.1-11     2023-04-25 [1] CRAN (R 4.2.0)
##  forcats      * 1.0.0      2023-01-29 [1] CRAN (R 4.2.0)
##  foreach        1.5.2      2022-02-02 [1] CRAN (R 4.2.0)
##  fs             1.6.3      2023-07-20 [1] CRAN (R 4.2.0)
##  generics       0.1.3      2022-07-05 [1] CRAN (R 4.2.0)
##  ggplot2      * 3.4.3      2023-08-14 [1] CRAN (R 4.2.0)
##  ggsci          3.0.0      2023-03-08 [1] CRAN (R 4.2.0)
##  glue           1.6.2      2022-02-24 [1] CRAN (R 4.2.0)
##  gtable         0.3.4      2023-08-21 [1] CRAN (R 4.2.0)
##  highr          0.10       2022-12-22 [1] CRAN (R 4.2.0)
##  hms            1.1.3      2023-03-21 [1] CRAN (R 4.2.0)
##  htmltools      0.5.6      2023-08-10 [1] CRAN (R 4.2.0)
##  htmlwidgets    1.6.2      2023-03-17 [1] CRAN (R 4.2.0)
##  httpuv         1.6.11     2023-05-11 [1] CRAN (R 4.2.2)
##  iterators      1.0.14     2022-02-05 [1] CRAN (R 4.2.0)
##  jquerylib      0.1.4      2021-04-26 [1] CRAN (R 4.2.0)
##  jsonlite       1.8.7      2023-06-29 [1] CRAN (R 4.2.0)
##  knitr          1.43       2023-05-25 [1] CRAN (R 4.2.2)
##  labeling       0.4.2      2020-10-20 [1] CRAN (R 4.2.0)
##  later          1.3.1      2023-05-02 [1] CRAN (R 4.2.2)
##  lattice        0.21-8     2023-04-05 [1] CRAN (R 4.2.0)
##  lifecycle      1.0.3      2022-10-07 [1] CRAN (R 4.2.0)
##  lubridate    * 1.9.2      2023-02-10 [1] CRAN (R 4.2.0)
##  magrittr       2.0.3      2022-03-30 [1] CRAN (R 4.2.0)
##  MASS           7.3-60     2023-05-04 [1] CRAN (R 4.2.2)
##  Matrix         1.5-4.1    2023-05-18 [1] CRAN (R 4.2.0)
##  memoise        2.0.1      2021-11-26 [1] CRAN (R 4.2.0)
##  mime           0.12       2021-09-28 [1] CRAN (R 4.2.0)
##  miniUI         0.1.1.1    2018-05-18 [1] CRAN (R 4.2.0)
##  munsell        0.5.0      2018-06-12 [1] CRAN (R 4.2.0)
##  myTAI        * 1.0.1.9000 2023-08-29 [1] Github (drostlab/myTAI@699b78f)
##  pillar         1.9.0      2023-03-22 [1] CRAN (R 4.2.0)
##  pkgbuild       1.4.2      2023-06-26 [1] CRAN (R 4.2.0)
##  pkgconfig      2.0.3      2019-09-22 [1] CRAN (R 4.2.0)
##  pkgload        1.3.2.1    2023-07-08 [1] CRAN (R 4.2.0)
##  prettyunits    1.1.1      2020-01-24 [1] CRAN (R 4.2.0)
##  processx       3.8.2      2023-06-30 [1] CRAN (R 4.2.0)
##  profvis        0.3.8      2023-05-02 [1] CRAN (R 4.2.0)
##  promises       1.2.1      2023-08-10 [1] CRAN (R 4.2.2)
##  ps             1.7.5      2023-04-18 [1] CRAN (R 4.2.0)
##  purrr        * 1.0.2      2023-08-10 [1] CRAN (R 4.2.2)
##  R6             2.5.1      2021-08-19 [1] CRAN (R 4.2.0)
##  Rcpp           1.0.11     2023-07-06 [1] CRAN (R 4.2.0)
##  readr        * 2.1.4      2023-02-10 [1] CRAN (R 4.2.0)
##  remotes        2.4.2.1    2023-07-18 [1] CRAN (R 4.2.2)
##  rlang          1.1.1      2023-04-28 [1] CRAN (R 4.2.0)
##  rmarkdown      2.24       2023-08-14 [1] CRAN (R 4.2.0)
##  rstudioapi     0.15.0     2023-07-07 [1] CRAN (R 4.2.0)
##  sass           0.4.7      2023-07-15 [1] CRAN (R 4.2.0)
##  scales         1.2.1      2022-08-20 [1] CRAN (R 4.2.0)
##  sessioninfo    1.2.2      2021-12-06 [1] CRAN (R 4.2.0)
##  shiny          1.7.5      2023-08-12 [1] CRAN (R 4.2.0)
##  stringi        1.7.12     2023-01-11 [1] CRAN (R 4.2.0)
##  stringr      * 1.5.0      2022-12-02 [1] CRAN (R 4.2.0)
##  survival       3.5-7      2023-08-14 [1] CRAN (R 4.2.0)
##  tibble       * 3.2.1      2023-03-20 [1] CRAN (R 4.2.0)
##  tidyr        * 1.3.0      2023-01-24 [1] CRAN (R 4.2.0)
##  tidyselect     1.2.0      2022-10-10 [1] CRAN (R 4.2.0)
##  tidyverse    * 2.0.0      2023-02-22 [1] CRAN (R 4.2.0)
##  timechange     0.2.0      2023-01-11 [1] CRAN (R 4.2.0)
##  tzdb           0.4.0      2023-05-12 [1] CRAN (R 4.2.2)
##  urlchecker     1.0.1      2021-11-30 [1] CRAN (R 4.2.0)
##  usethis        2.2.2      2023-07-06 [1] CRAN (R 4.2.0)
##  utf8           1.2.3      2023-01-31 [1] CRAN (R 4.2.0)
##  vctrs          0.6.3      2023-06-14 [1] CRAN (R 4.2.0)
##  vroom          1.6.3      2023-04-28 [1] CRAN (R 4.2.0)
##  withr          2.5.0      2022-03-03 [1] CRAN (R 4.2.0)
##  xfun           0.40       2023-08-09 [1] CRAN (R 4.2.2)
##  xtable         1.8-4      2019-04-21 [1] CRAN (R 4.2.0)
##  yaml           2.3.7      2023-01-23 [1] CRAN (R 4.2.0)
## 
##  [1] /Library/Frameworks/R.framework/Versions/4.2/Resources/library
## 
## ──────────────────────────────────────────────────────────────────────────────